Simulations of nonhelical hydromagnetic turbulence 
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Nonlielical hydromagnetic forced turbulence is investigated using large scale simulations on up 
to 256 processors and 1024'' meshpoints. The magnetic Prandtl number is varied between 1/8 
and 30, although in most cases it is unity. When the magnetic Reynolds number is based on 
the inverse forcing wavenumber, the critical value for dynamo action is shown to be around 35 
for magnetic Prandtl number of unity. For small magnetic Prandtl numbers we find the critical 
magnetic Reynolds number to increase with decreasing magnetic Prandtl number. The Kazantsev 
k^^'^ spectrum for magnetic energy is confirmed for the kinematic regime, i.e. when nonlinear effects 
are still unimportant and when the magnetic Prandtl number is unity. In the nonlinear regime, the 
energy budget converges for large Reynolds numbers (around 1000) such that for our parameters 
about 70% is in kinetic energy and about 30% is in magnetic energy. The energy dissipation rates are 
converged to 30% viscous dissipation and 70% resistive dissipation. Second order structure functions 
of the Elsasser variables give evidence for a k~^^^ spectrum. Nevertheless, the three-dimensional 
spectrum is close to k~^^'^ , but we argue that this is due to the bottleneck effect. The bottleneck 
effect is shown to be equally strong both for magnetic and nonmagnetic turbulence, but it is far 
weaker in one-dimensional spectra that are normally studied in laboratory turbulence. Structure 
function exponents for other orders are well described by the She-Leveque formula, but the velocity 
field is significantly less intermittent and the magnetic field is more intermittent than the Elsasser 
variables. 



I. INTRODUCTION 



Dynamo action, i.e. the conversion of kinetic energy 
into magnetic energy, plays an important role in astro- 
physical bodies ranging from stars to galaxies and even 
clusters of galaxies. The gas in these bodies is turbu- 
lent and the magnetic Reynolds numbers are huge (10^° 
to 10^"). This suggests that there should be dynamo 
action and that the magnetic fields should be amplified 
on the time scale of the turbulent turnover time [3, Q . 
The magnetic fields of many astrophysical bodies show a 
great deal of spatio-temporal order (e.g. an 11-year cy- 
cle and equatorward migration of the magnetic field in 
the case of the sun), which can basically be explained by 
the helicity effect |^ H, Q. In many other astrophysi- 
cal environments, however, the helicity effect is probably 
completely irrelevant (e.g. in the solar wind 0, clusters 
of galaxies 7] , and the early universe after recombination 
0). In all these cases the magnetic field does not show 
spatio-temporal order of the type known for helical hy- 
dromagnetic turbulence Both analytic theory ,9| as 
well as simulations 10] have long shown that dynamo ac- 
tion is possible even without kinetic helicity, but that the 
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field is then spatially highly intermittent with substan- 
tial power at small length scales. Hydromagnetic tur- 
bulence has recently also been studied in the laboratory 
[Til IT^j ] , but there the magnetic Reynolds numbers 
are still rather small. Therefore, numerical simulations 
are currently the most powerful tool. 

In spite of significant progress over the past two 
decades, the form of the energy spectrum at large mag- 
netic Reynolds numbers is still a matter of debate. Par- 
ticular progress has been made in the case where there 
is a large scale field. Goldreich & Sridhar 14] have pro- 
posed that the magnetic energy spectrum has a k~^^^ 
inertial range due to the anisotropy imposed by the lo- 
cal magnetic field. This was in conflict with the ear- 
lier Iroshnikov-Kraichnan k~^/'^ spectrum, which as- 
sumes isotropy, but can now be ruled out 0, 0, 0] . 
Some simulations still show a k~^/'^ spectrum, but this is 
probably due to the bottleneck effect; see below. A k~^^^ 
total energy spectrum has however clearly been seen in 
decaying hydromagnetic turbulence Here the total 
energy spectrum is simply the sum of kinetic and mag- 
netic energy spectra. 

The case of an imposed large scale magnetic field with 
approximate equipartition strength is in some respects 
similar to the case of a self-generated large scale field 
which emerges if there is helicity in the flow [J, 0, • In 
the latter case of helical turbulence, and for unit magnetic 
Prandtl number, kinetic and magnetic energy spectra are 
in almost perfect equipartition on all scales smaller than 
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the forcing scale (see Fig. 11 of Ref. |^), and the two spec- 
tra tend to approach a k"^/'^ inertial range as the mag- 
netic Reynolds number is increased. In the former case, 
equipartition throughout the inertial range may require 
that the energy of the imposed field is comparable to the 
kinetic energy; for stronger fields the magnetic energy 
becomes suppressed at small scales If the magnetic 
Prandtl number is larger than unity, the resistive cutoff 
is p rolonged to larger wavenumbers by a spectrum 
plf . Nevertheless some authors have suggested that even 
in the nonlinear regime the spectral magnetic energy in- 
creases toward smaller scales, similar to the kinematic 
regime where the energy spectrum scales as fc"*"^/^ [^12^. 

In the absence of an imposed large scale magnetic field, 
and with no helicity, the situation is in many ways differ- 
ent. Early simulations suggest that the magnetic field is 
dominated by small scale power ilOj. Even for a magnetic 
Prandtl number of unity the magnetic energy exceeds the 
kinetic energy at small scales 2'^, 25] , a result that is 
otherwise (with imposed field or with helicity) only ob- 
tained for magnetic Prandtl numbers larger than unity 

The main problem in determining the energy spectrum 
at large hydrodynamic and magnetic Reynolds numbers 
is the lack of a proper inertial range in the magnetic 
field. Another problem is that in the absence of helicity 
the dynamo is much weaker and one needs significantly 
(about 20 x) larger magnetic Reynolds numbers before 
the dynamo is even excited. In practice, this means that 
a resolution of 128'^ or more is mandatory. 

The shortness or even lack of an inertial range has 
frequently led to the use of hyperviscosity and hyperdif- 
fusivity (fol |2^ , which has the tendency to extend the 
inertial range and to shorten the dissipative subrange. 
However, in recent years it has become clear that these 
modifications to the viscosity and diffusion operators can 
affect major parts of the inertial subrange and make it 
shallower. This is also referred to as the bottleneck ef- 
fect which is present already for ordinary viscosity |27j . 
but it becomes greatly exaggerated with hyperviscosity; 
compare Figs 8 and 11 of Ref. helical dynamos 

the saturation time becomes significantly prolonged and 
the final saturation level is artificially enhanced 28] . It is 
therefore important to use simulations with regular vis- 
cosity and magnetic diffusivity at high enough resolution. 

A lot of work has already been done in order to deter- 
mine the form of the magnetic energy spectrum at large 
magnetic Reynolds numbers. Some of the relevant papers 
are listed in Table where we indicate the main prop- 
erties of their set-ups. In addition to the papers men- 
tioned above we have included a number of additional 
ones m m m m m Hi Hi, some of which will be 
discussed below in more detail. 

In a recent paper we have already presented some ini- 
tial results on energy spectra at a resolution of up to 
1024'^ meshpoints; see Ref. [23> hereafter referred to as 
Paper I. In the present paper we discuss the associated 
results for the critical magnetic Reynolds number for dy- 



TABLE I: Summary of earlier work on hydromagnetic turbu- 
lence simulations. In addition to the reference of each paper 
we give for convenience also a more descriptive abbreviation. 
In Ref. the largest resolution was A''^ — 256^ x 512, so 
we have listed the geometrical mean corresponding to 322'^. 
The column 'hyper' indicates whether or not hyperviscosity 
or hyperdiffusivity have been used; 'artif stands for shock- 
capturing artificial viscosity. In the column 'forced', the ab- 
breviation 'hel' stands for helical, 'noh' for non-helical, 'decay' 
for decaying turbulence without forcing, and '~noh' means 
that there is some small fraction of helicity, but due to a small 
scale separation between the forcing and the box scale the 
simulation is practically non-helical. In the column (iJ), 'ext' 
indicates the use of an external field. The column 'compr' 
signifies whether the models are compressible. 



Paper 


Ref. 


N hyper 


forced 


{B) 


compr 


PrM 


MFP81 


m 


64 


both 


hel/noh 





no 


1 


BOX 




120 


no 


hel 





yes 


0.1-100 


KYM91 


[23j 


128 


yes 


hel 





no 


1 


MMDMOl 


[29] 


128 


no 


decay 


ext/0 


no 


1 


OMG98 


[30] 


128 


no 


decay 


ext/0 


both 


1 


CL03 


[31J 


216 




decay 


ext 


yes 


1 


MLKB98 


[32J 


256 


artif 


decay 





yes 


1 


CVOO 


m 


256 


both 


~noh 


ext 


no 


1 


MCOl 


[25j 


256 


both 


noh 


ext/0 


no 


1-2500 


MGOl 


[18j 


322 


yes 


noh 


ext 


no 


1 


CLV02 


[33J 


384 


yes 


~noh 


ext 


no 


large 


BNP02 


[34] 


500 


yes 


noh 


ext 


yes 


1 


BMOO 


m 


512 


both 


decay 





no 


1 


HBD03 




1024 


no 


noh 





yes 


1 



namo action, we consider a range of different values of the 
magnetic Prandtl number, compare with the case of finite 
magnetic helicity, and present visualizations of the mag- 
netic field and the dissipative structures in hydromag- 
netic turbulence. We look in detail at the energy spectra 
and the structure functions of total, kinetic and magnetic 
energies. The simulations have been carried out using 
the Pencil Code [s^ which is a memory-efficient high- 
order finite difference code using the 2A^-RK3 scheme of 
Williamson 0]. Our approach is technically similar to 
that used in Ref. 0]. The Pencil Code is fully compress- 
ible; we therefore consider the weakly compressible case 
(the Mach number is around 0.1), which can be consid- 
ered as an approximation to the incompressible case. A 
list of parameters for most runs discussed in the present 
paper is given in Table Ull 



II. EQUATIONS 

We adopt an isothermal equation of state with con- 
stant (isothermal) sound speed Cg, so the pressure p is 
related to the density p hy p = pc^. The equation of 
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TABLE II; Summary of runs with PrM = 1 (thus Re — Rbm) 
and forcing at k{ = 1.5. In all cases, except D2, we have 
Mrms ~ 0.12. Run D2 had a stronger forcing, and there- 
fore Urms ~ 0.18. The Taylor microscale Reynolds number is 
Rca = uid£\/i-', where lilD - — 'firms /VS is the rms velocity in 
one direction and £x = V^Wrms/i^rms is the Taylor microscale, 
where uj^-ms is the rms vorticity. 



Run Resolution 


u X 10* 


Rca 


Rcm 


^rms 


eK X 10"* 


EM X lO"* 


A 64^ 


7.0 


80 


120 


0.052 


1.0 


1.2 


B 128^ 


4.0 


110 


190 


0.060 


0.79 


1.4 


C 256^ 


2.0 


160 


420 


0.062 


0.78 


1.6 


D 512^ 


1.5 


190 


540 


0.072 


0.68 


1.7 


D2 512^ 


2.0 


180 


600 


0.092 


2.2 


5.0 


E 1024^ 


0.8 


230 


960 


0.075 


0.63 


1.5 



motion is written in the form 

J xB 



Dm 2^, 



(1) 



where D/Di = d/dt + m • V is the advective derivative, 
= V X JB//io is the current density, /io is the vacuum 
permeability, 

Fvisc = (V^tt + iVV • M + 2S • Vlnp) (2) 

is the viscous force, u — const is the kinematic viscosity, 



1 / dui_ 

2 V 9a;,- 



dxi 



V • u 



(3) 



time-integrated forcing function to be independent of the 
length of the time step 5t, the normalization factor N has 
to be proportional to St'^^"^. On dimensional grounds it 
is chosen to be = /oCs(|fc|cs/(5i)^/^, where /o is a nondi- 
mensional forcing amplitude. The value of the coefficient 
/o is chosen such that the maximum Mach number stays 
below about 0.5; in practice this means /o = 0.02 . . . 0.05. 

At each timestep we select randomly one of many pos- 
sible wavevectors in a certain range around a given forc- 
ing wavenumber. The average wavenumber is referred to 
as fcf. We force the system with nonhelical transversal 
waves, 



fk = {kx e) /\/fc2 - (fc . e)2 



(7) 



where e is an arbitrary unit vector not aligned with fe; 
note that |/fcp — 1. 

The resulting flows are characterized by the kinetic and 
magnetic Reynolds numbers, 



Re 



Re 



M 



(8) 



respectively. Their ratio is the magnetic Prandtl number, 
PrM = = ReM/Re, (9) 

which is unity for most of the runs. 

We use non-dimensional quantities by measuring 
length in units of l/ki (where ki = 2tt/L is the smallest 
wavenumber in the box of size L), speed in units of the 
isothermal sound speed Cg, density in units of the initial 
value po, and magnetic field in units of (/iopoCs)^^^- 



is the traceless rate of strain tensor, and / is a random 
forcing function that consists of non-helical plane waves 
(see below). The continuity equation is written in terms 
of the logarithmic density, 



Dlnp 



(4) 



and the induction equation is solved in terms of the mag- 
netic vector potential A, with B — 'V x A, so 



dA 

'dt 



ux B + r/V^A, 



(5) 



where rj = const is the magnetic diffusivity. 

We use periodic boundary conditions in all three di- 
rections for all variables. This implies that the mass in 
the box is conserved, i.e. (p) = pq, where po is the value 
of the initially uniform density, and angular brackets de- 
note volume averages. We adopt a forcing function / of 
the form 

fix, t) = Yie{Nfk(t) exp[ifc(i) • x + i0(t)]}, (6) 

where x is the position vector. The wave vector k{t) 
and the random phase — tt < (/)(i) < tt change at ev- 
ery time step, so f{x,t) is 5-correlated in time. For the 



III. THE KINEMATIC PHASE AND 
APPROACH TO SATURATION 

As initial condition we use a weak random seed mag- 
netic field. In this section we consider the time interval 
during which the amplitude of the magnetic energy spec- 
trum grows exponentially and is still small compared to 
the kinetic energy spectrum for all wavenumbers. 



A. Critical magnetic Reynolds number 

We have determined the growth rate of the rms mag- 
netic field, A = dlni3rms/dt, for different values of Rcm 
and PrM. Interpolating the curves A = A(ReM) through 
zero, we find the critical value, Re^'^'*\ which is roughly 
independent of fct; see Fig. ^ For the case PrM = 1 



we find Re 



(crit) 
M 



35, which is consistent with the value 



obtained using a modified version of the eddy-damped 
quasi-normal Markovian (EDQNM) approximation [Sg] 

which gives Rel^"*^ = 29. 

As the value of PrM is lowered, the critical magnetic 
Reynolds number increases slowly like 



Re 



(crit) 
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35 Pr 



-1/2 



M 



(10) 
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FIG. 1: Growth rate versus ReM for different values of kf, 
and for PrM = 1- The curves represent linear fits through 
the data. Note that the critical value is around 35 for all the 
different runs. The resolution varies between 64^ and 256'', 
and /o = 0.05 in all runs, resulting in Urms ~ 0.2. 



sists to very small values of PrM that are relevant to liq- 
uid metal experiments or to stellar convection zones. We 
note, however, that the EDQNM approximation predicts 
that the critical magnetic Reynolds number is inde- 
pendent of PrM for a large range of Re and PrM ■ Based on 
application of the Kazantsev model, Schekochihin et al. 
|39l | have shown that, if the correlation time is assumed to 
be independent of wavenumber, there exists a finite value 
of PrM below which dynamo action is impossible. On the 
other hand, if the correlation time is proportional to the 
eddy turnover time, ~ k~^^^, dynamo action should be 
possible even for very small values of PrM ■ 

For larger magnetic Reynolds numbers the growth rate 
of the magnetic field approaches a Re]^^ dependence; see 
Fig. 131 Such a power law would be expected if the growth 
rate is proportional to the eddy turnover time at the 
dissipation wavenumber, fed, i.e. A oc fc^'^ cx Re]^^, where 

we have used Rcm oc k^^^ We must emphasize, 

however, that the Rcm scaling is as yet rather short. 
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FIG. 2: Critical magnetic Reynolds number as a function of 
magnetic Prandtl number for runs with different forcing scale. 
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FIG. 3: Log- log plot of growth rate versus ReM for fcf — 1.5, 
and PrM ~ 1- Note that A/(«rmsfcf) seems to approach 

1/2 

Rcj^j /300 for large values of Rcm. The resolution varies be- 
tween 128^ and 512=*. 



see Fig. |21 It remains uncertain whether this scaling per- 



B. The Kazantsev spectrum 

Under the somewhat unrealistic assumption that the 
velocity field is i5-correlatcd in time, one can, for a given 
spatial correlation function of the velocity, and ignor- 
ing magnetic feedback, derive an evolution equation for 
the correlation function of the magnetic field 0, |4l| (see 
Ref. ^3 for a nonlinear extension). This can be rewrit- 
ten as an integro-differential equation in fc-space which, 
in turn, can be written as a diffusion equation in fc-space 
if the velocity field has only power at large scales j3, 1231 ■ 
The result is E^{k) cx k^/^K„(^x)ik/k,^), where is 
the Macdonald function (modified Bessel function of the 
third ki nd) , A ~ 3/4 is an eigenvalue, and n(X) « 0; 
see Ref. [i^ for details. In linear theory the amplitude 
of the solution is however undetermined and grows ex- 
ponentially if the magnetic Reynolds number exceeds a 
critical value of about 30-60 [ij ■ 

In comparison with our simulations we see qualitative 
agreement as far as the A:'^/^ slope is concerned. As long 
as the magnetic field is weak, the velocity field has an 
energy spectrum with the expected k~^/^ Kolmogorov 
scaling; see Fig0] During this phase the spectral mag- 
netic energy grows at all wavenumbers exponentially in 
time and the spectrum has the expected k^^^ Kazantsev 
slope; see Fig0] The convergence toward Kazantsev 
scaling is seen in Fig. |31 where we have plotted the mag- 
netic energy spectrum for runs with Rcm between 120 
and 540. (The kinematic growth phase for Run E was 
not available, because we restarted Run E from Run D 
after it had already reached saturation.) 

Originally, Kazantsev obtained the k^^^ spectrum un- 
der the assumption that the velocity has power only at 
large scales, which would correspond to a large value of 
PrM- Simulations for large PrM have indeed confirmed 
the Kazantsev slope . Our results now show that the 
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FIG. 4: Early spectra of kinetic and magnetic energy, nor- 
malized by iu^ms/'^i, during the kinematic stage of run D2. 
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FIG. 5: Convergence of the magnetic power spectrum toward 
the k^^^ scaling as ReM increases. All spectra correspond to 
the kinematic stages of Runs A-D with Prm = 1, but have 
been rescaled to make them coincide at small values of k. 
The dimension on the ordinate is therefore arbitrary. 



Kazantsev spectrum is also obtained for PrM — 1. 

The values of 77 (= v) and ReM used for the different 
runs discussed above are summarized in Table HTl 



C. Approach to saturation 

One would naively expect that the onset of satura- 
tion happens rapidly on a dynamical time scale. That 
this is not the case in helical hydromagnetic turbulence 
came originally as a surprise, but the reason is now well 
understood to be a consequence of asymptotic magnetic 
helicity conservation. The same argument does not ap- 
ply to nonhelical turbulence. Nevertheless, there may 
still be a slow-down in the saturation behavior at small 
wavenumbers 

The clarification of this question is hampered by the 
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FIG. 6: Saturation behavior of the spectral magnetic energy 
at wavenumbers fc = 1 (solid line) and fc = 16 (dashed line) 
for run D2. Note the slow saturation behavior for k = 1. 



fact that one needs very high resolution before one can 
with certainty distinguish resistive time scales from dy- 
namical ones. Figure IHI shows the saturation behavior 
for two different values of fc in a run with 512'^ mesh- 
points. The initial growth is clearly exponential both 
at /c = 1 and at fc = 16. However, when the magnetic 
energy at fc = 16 saturates at i ~ 100/ (wrmsfcf ), the mag- 
netic energy at fc = 1 seems to continue growing approx- 
imately linearly "S]; see Fig. for 100 < Mjmsfcfi < 200. 
Looking at the simulation at later times one sees, how- 
ever, that the time sequence is actually very bursty 
and that the approximately linear growth in the inter- 
val 100 < Urmsfcfi < 200 was only a transient. (This is 
best seen in the inset of Fig. |B| where the energies are 
shown in a linear plot.) We can therefore not confirm 
with certainty that the difference in the saturation times 
for fc = 1 and fc = 16 is explained by the difference be- 
tween dynamical and resistive time scales. 

In the following section we present further properties 
of these runs after the time when the field has reached 
saturation. 



IV. DYNAMICALLY SATURATED PHASE 

When the magnetic energy has reached a certain frac- 
tion of the kinetic energy, the magnetic field stops grow- 
ing exponentially and eventually reaches a statistically 
steady state. In this section we discuss the properties of 
this state. 



A. Energy spectra 

As saturation sets in, the spectral magnetic energy be- 
gins to exceed the spectral kinetic energy at small scales; 
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FIG. 7: Three-dimensional (upper panel) and one- 
dimensional (lower panel) time averaged spectra, normalized 
by |u?ms/^ii for Run E with a resolution of 1024'^ meshpoints. 
The energy spectra have been averaged over a period of five 
turnover times (urmsfcf)"^- 



see the upper panel of Fig. [7| where we also show the 
total energy spectrum, = + E^. There is a short 
inertial range with an approximate Kolmogorov k~^^'^ 
spectrum, but, as already reported in Paper I and dis- 
cussed in Ref. 0|j there is a strong bottleneck effect 
in the three-dimensional spectra which is less strong in 
one-dimensional spectra. Therefore we plot in Fig. [T] for 
comparison both three-dimensional and one-dimensional 
spectra for the same Run E. With one-dimensional spec- 
tra we mean spectra calculated from variations along one 
of the coordinate directions only. Throughout this paper, 
we refer to the one-dimensional spectra as the sum of the 
longitudinal plus two times the transversal spectra p^ . 
To obtain the three-dimensional spectra we have inte- 
grated over three-dimensional shells in k space, assuming 
isotropy. 

Much of the large-scale kinetic energy is probably 
transferred directly to smaller-scale magnetic fields. At 
large scales the three-dimensional magnetic energy spec- 
trum is weakly increasing, it peaks at fc = 5, and then 
joins the k~^^^ slope of kinetic energy, but with an 



approximately 2.5 times larger amplitude. The one- 
dimensional energy spectrum, on the other hand, is 
monotonically decreasing also for small wavenumbers. 

All spectra terminate around the nominal dissipation 
cutoff wavenumber, = {e/v^Y^'^, where v = rj and 
e = eK + Cm is the total energy dissipation rate per unit 
mass, and 

eK = 2j/(pS^)/po and cm = Wo(J'^)/po (H) 

are the contributions from viscous and ohmic dissipation, 
respectively. As the resolution is increased, the inertial 
range of the total energy spectrum becomes progressively 
longer (see Paper I). For Run E we find k^ ~ 142. The 
'hook' in the spectrum for k > 400 (see Fig.[7|l is probably 
a consequence of finite resolution and is typical also of 
turbulence simulations using spectral codes j45|. 

For the numerical simulations to work and to repro- 
duce the turbulent cascade reliably, a certain minimum 
dynamic range in the compensated spectrum k^/^Ek{k) 
is necessary. We emphasize that this dynamic range in- 
creases with increasing Reynolds number, as can be seen 
in Fig. |S1 where we have used viscosities close to the min- 
imum values required for reliable results. 

Looking at one-dimensional spectra of the magnetic, 
kinetic and total energy (denoted by Sf''"", E^-'" 
and e'^'^^), compensated with e~'^^^k^^^, one sees that 
around k — 10 there is a short range where all three com- 
pensated spectra are flat; see Fig.O We emphasize again 
that this is seen only in the higher resolution runs. 

There is some ambiguity as to what one calls the in- 
ertial range. Here we refer to inertial range as only the 
range where the total energy spectrum is compatible with 
a k"^/^ slope. The kinetic energy alone, however, can 
show a k^^/^ slope already for smaller wavenumbers; see 
Fig.El 



B. Comparison with nonmagnetic turbulence 

For comparison we show the resulting spectra for 
purely hydrodynamic turbulence without magnetic fields 
and find again very similar bottleneck behavior, as shown 
in Fig. ^1 A similar bottleneck effect is also found in 
other numerical simulations on up to 1024^ meshpoints 
[4^ lirj . As one increases the resolution even further 
(up to 4096'^ meshpoints) the bottleneck assumes an 
asymptotic shape, and begins to separate from the in- 
ertial range. A weak bottleneck effect is found even in 
wind tunnel experiments |49l ISOf . In such experiments 
one usually measures one-dimensional longitudinal en- 
ergy spectra, and hence a much weaker bottleneck effect 
is expected |4J1 . A physical explanation of the bottleneck 
effect is found in Ref. [27j . 
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FIG. 8: Convergence of compensated three-dimensional (top) 
and one- dimensional (bottom) energy spectra of kinetic en- 
ergy. The vertical arrows at the top show the dissipation 
cutoff wavenumber {e/u^y^'^ for the different resolutions. 



C. Subinertial range behavior 

For the saturated state, Fig. ^2 shows that for A: < 5 
the magnetic energy spectrum follows approximately a 
k^^^ behavior, as was originally expected by Batchelor 
0. The same slope has also been found for convective 
turbulence at the time when the magnetic field is still 
weak 51], and, more recently (b ut p robably for different 
reasons) for ABC- flow dynamos l53|. 

For wavenumbers much less than the forcing wavenum- 
ber (subinertial range), the magnetic and kinetic energy 
spectra increase with k approximately like k^. The the- 
ory for this spectrum is reviewed in the book by Lesieur 
|5a| . This scaling is best seen in spectra where the tur- 
bulence is forced a,t kf ^ ki = 1; see Fig. 1121 where 
kf = 15. In this figure the kinetic spectrum is somewhat 
shallower and scales more like k^'^. In this run, the mag- 
netic energy is saturated, and yet it is only about 4% 
of the kinetic energy. This is mainly because the mag- 
netic Reynolds number is only about 37, so the run is 
just weakly supercritical. 
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FIG. 9: Compensated magnetic, kinetic and total three- and 
one- dimensional energy spectra for Run E. 



D. Reynolds number dependence 

In Fig. El we show total energy spectra for three val- 
ues of Re = ReM between 270 and 960. In the first 
case with Re = Rcm = 270 the compensated spectrum 
shows a reasonably flat range with a Kolmogorov con- 
stant Ckym ~ 1-3, which is somewhat less than the value 
of 2.1 found by Kida et al. @|. In the second case with 
Re = Rcm = 440 the compensated spectrum shows clear 
indications of excess power just before turning into the 
dissipative subrange. This is just the bottleneck effect 
and it becomes even more dramatic in the third case 
with Re = Rcm = 960. Thus, the bottleneck effect is 
only seen at sufficiently large resolution. 

In order to assess the reliability of the results we have 
carried out a convergence study using the same value of 
v, but different mesh resolution; see the second panel 
of Fig. El In both cases, v = r] = 2 x lO"'*, while 
M„„s = 0.13 for 512^ meshpoints and 0.12 for 256^. The 
energy dissipation is also similar, e = 2.8 x lO""' for 512'^ 
meshpoints and 2.3 x 10"^ for 256^. The compensated 
energy spectra agree quite well for the two different res- 
olutions, and both show excess power (bottleneck effect) 
just before turning into the dissipative subrange. This 
supports the conclusion that the excess power is not a 
numerical artifact. 
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k/kx 

FIG. 10: Same as Fig.m but for a run without masnetic fields 
(pure hydrodynamics), for three different values of Re. Note 
the appearance of a bottleneck effect in the range 10 < fc < 50, 
particularly in the 3-D spectra. For these runs fcf = 1.5, and 

Itrms ~ 0.18. 
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FIG. 11: Saturated three-dimensional magnetic energy spec- 
tra compensated by k~^^^ . The results are comparable with 
the vorticity-like k^^^ scaling at small wavenumbers (fc < 4). 
The vertical arrows at the top show the dissipation cutoff 
wavenumber {e/v'^Y^'^ for the different resolutions. 
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FIG. 12: The sub-inertial range of hydromagnetic turbulence 
forced at fcf = 15 and normalized by ^u^ui^/ki. Rcm ~ 37, 
256^ meshpoints, PrM ~ 1. 



E. Convergence of energy and dissipation rate 

In the saturated state, the fractional magnetic and ki- 
netic energies tend to a constant value at large Reynolds 
number, with 

Em ■■ Ek « 0.3 : 0.7, (12) 

so Em/Ek is about 0.4; see the upper panel of Fig. 1141 
This fraction may still depend on the forcing wavenumber 
since the infrared part of the kinetic and magnetic energy 
spectra can have different slopes; see Fig. ^1 Naively 
one would just compare Em/Ek for the three non-helical 
simulations in Fig. and find the ratio to decrease as 
the forcing wavenumber is increased. This is however 
not correct since all these runs have different (magnetic) 
Reynolds numbers. A more comprehensive study, where 
the Reynolds number is kept constant, is therefore neces- 
sary in order to find whether or not Em/E^ really does 
depend on the forcing wavenumber. 

Equally important is the fact that also the energy dis- 
sipation rates are converged for large Reynolds numbers, 
such that 

CM : EK » 0.7 : 0.3, (13) 

i.e. Em exceeds ck by a factor of about 2.3; see the two 
lower panels of Fig. E| The reciprocal correspondence 
of the ratios in Eqs (|12|) and H13|) is coincidental. 

The fact that the dissipation rates for both mag- 
netic and kinetic energies are asymptotically indepen- 
dent of Reynolds number is consistent with the basic 
Kolmogorov phenomenology that leads to the scale-free 
^-5/3 spectrum. This result seems to exclude the possi- 
bility that in the large Reynolds number limit the mag- 
netic energy spectrum peaks at small scales. It is worth 
mentioning that the ratio eyi/^K depends only weakly on 
PrM. In Table ITTll we show that a 15 fold increase of the 
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FIG. 13: Comparison of compensated three-dimensional to- 
tal energy spectra for runs with magnetic Prandtl number 
unity but with different Reynolds numbers. In all runs the 
horizontal dash-dotted hne represents the value 1.3. In the 
second panel two runs with the same Reynolds numbers, but 
different resolution are compared. 



TABLE III: Energy dissipation rates for four runs with dif- 
ferent magnetic Prandtl numbers PrM, showing that eu/^K is 
only weakly dependent on PrM- 
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magnetic Prandtl number decreases the energy dissipa- 
tion ratio only by a factor of about two; the data can be 
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parameterized by the power law em/^K ~ 2.2 Pr 

About 70% of the energy that goes into the Kol- 
mogorov cascade is eventually converted into magnetic 
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FIG. 14: Relative magnetic and kinetic energy (upper panel), 
their respective relative dissipation rates (middle panel), and 
the energy dissipation rates in units of eo = kiu'}^^ (last 
panel), as a function of Reynolds number for PrM = 1- 



energy by dynamo action and is then finally converted 
into heat via resistivity. A sketch of the energy budget 
is given in Fig. 1151 

In the present simulations the thermal energy bath is 
not explicitly included, but we note that in hydromag- 
netic turbulence simulations with shear and rotation the 
resistive heat can become so important that the temper- 
ature can increase by a factor of 10; see Ref. [s^ . 



F. Large magnetic Prandtl numbers 

It has previously been argued that for PrM ^ 1 
the magnetic energy spectrum is peaked at small scales 
26, 55, 56J. For PrM — 1 this claim cannot be con- 
firmed at large Reynolds number (see Paper I). The orig- 
inal motivation for a peak at small scales is based on lin- 
ear theory 9], which predicts a /c"*"^/^ spectrum; see also 
Sec. lIIIBI However, the original Kazantsev model is only 
valid in the limit where the velocity has only large scale 
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FIG. 15: Schematic view of the energy transport. Most of 
the energy (70%) resides in the kinetic energy reservoir, but 
only 30% of the total energy input is dissipated directly by 
viscous heating. Instead, 70% of the energy flows into the 
magnetic energy reservoir, and is finally dissipated by Ohmic 
heating. 



FIG. 17: Autocorrelation functions of magnetic field and 
velocity. Note that the autocorrelation functions are nearly 
independent of resolution and Reynolds number. The velocity 
correlation length is ~ 3 while the magnetic correlation length 
is ~ 0.5. 



^ 10-2 




FIG. 16: Magnetic energy spectra for runs with magnetic 
Prandtl numbers ranging from 0.3 to 30. 



components, which corresponds to Pim ^ 1. In order to 
see how our results change with varying magnetic Prandtl 
number we have calculated models for different values of 
PrM- One of the results is the possible emergence of a 
tail in the magnetic energy spectrum; see Fig. 1161 
The k^^ tail has recently been found in large PrM sim- 
ulations with an imposed magnetic field [2l|. The k^^ 
spectrum has its roots in early work by Batchelor jS?] for 
a passive scalar and MofFatt [Sg for the magnetic case. 

In the run with PrM — 30 the viscous cutoff wavenum- 
ber is ki, = {ck/v^Y^'^ ~ 12, so the k~^ tail is expected 
for wavenumbers larger than that. The plot seems to 
suggest that the entire inertial range could have a 
spectrum, although this may well be an artifact of an in- 
sufhcient inertial range. Instead, a more likely scenario is 
that for large hydrodynamic Reynolds numbers and large 



magnetic Prandtl numbers there is still a range for 

both kinetic and magnetic energies, followed by a 
subrange for magnetic energy beyond the viscous cutoff 
wavenumber. In any case, the peak of magnetic energy 
would still be at small wavenumbers. In summary, there- 
fore, we find no indication of a peak of the magnetic 
energy spectrum at the resistive wavenumber. 



G. Autocorrelation functions 



In the context of the Zeldovich |59|, |60| stretch-twist- 
fold dynamo, the shape of the auto-correlation function 
of the magnetic field, 



WB{r) = {B{x)-B{x+r))/{B^ 



(14) 



plays an important role; see also Ref. [61|. We have 
calculated wb(^) from the Fourier transform of the 
three-dimensional, time-averaged magnetic energy spec- 
tra, E^{k); see Fig. El The diffusive scale corresponds 
to the thickness of the narrow spike of ^^(r) around the 
origin. The typical scale over which the magnetic field 
changes direction is the correlation length which corre- 
sponds to the scale where WBir) has its minimum. We 
have compared WB{r) with the velocity autocorrelation 
function, Wuir), which is defined in an analogous man- 
ner. Note that, because of isotropy, Wuir) and WB{r) are 
only functions of r = \r\. Similar autocorrelation hmc- 
tions have also been seen in simulations of convective 
dynamos (sj . Figure El shows that the velocity corre- 
lation length is 3 while the magnetic field correlation 
length is ~ 0.5. Clearly, the magnetic correlation length 
is much shorter than the velocity correlation length, but 
it is practically independent of Re (= Rcm) and certainly 
much longer than the resistive scale, ~ 27r/fcd ~ 0.04. 
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H. Structure functions 

In numerical turbulence the signed structure functions 
for odd moments are usually not well converged. It is 
therefore customary to use unsigned structure functions, 
defined as 

S,il) = {\zix+l)-zixr), (15) 

where z{x) is one of the two Elsasser variables 

z^^u±B/^/^. (16) 

The structure function exponents Sp{l) normally show a 
power-law scaling 



Sp{l) oc f^, 



(17) 



where (^p is the pth order structure function scaling ex- 
ponent. Here, because of isotropy, we usually consider 
the dependence on one spatial coordinate, x, but we also 
discuss the more general case below. 

In Paper I we have shown that ^3 « 1.0 and ^4 w 1.3. 
This is our strongest evidence that the asymptotic in- 
ertial range scahng is and that the k~^^^ scaling 

seen in the upper panel of Fig.[7|is due to the bottleneck 
effect |44| . True Iroshnikov-Kraichnan scaling would im- 
ply C4 = 1, whereas = 1 is consistent with Goldreich- 
Sridhar scaling. 

We find that the second-order scaling exponent is (2 = 
0.7, which again indicates that the inertial range has a 
slope k~''^~^ ~ k~^/^, in support of our findings from 
the one-dimensional energy spectra. Making use of the 
extended self similarity hypothesis [G^I we find that for 
all p, the values of (p are consistent with the generalized 
She-Leveque formula |6^ 



p 

9 



C 



2/3\ 
C j 



p/3' 



(18) 



where C is interpreted as the codimension of the dissi- 
pative structures. We find that C — 2 (corresponding 
to 1-dimensional, tube-like dissipative structures) gives a 
reasonable fit to the longitudinal structure function ex- 
ponents; see Fig.^J If we allow for fractal codimensions, 
then C = 1.85 gives the best fit for the longitudinal struc- 
ture function exponents and C — 1.45 for the transversal 
ones. 

A similar difference between transversal and longitudi- 
nal structure functions has previously been found pials^ . 
It has been argued '63| that this difference is an artifact 
of the forcing being in the same direction as the direction 
in which the structure functions are calculated. In our 
simulation, however, the forcing is chosen randomly in 
an isotropic manner. In addition we have also performed 
a calculation were we average structure functions calcu- 
lated in 91 different directions distributed isotropically 
over the unit sphere. In this calculation we find again 
the same difference between longitudinal and transver- 
sal structure functions. Additional support to this result 
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FIG. 18: Longitudinal and transversal structure function ex- 
ponents for the Elsasser variables for run D. 



comes from the fact that longitudinal one-dimensional 
energy spectra are slightly steeper than the transversal 
ones. A possible explanation for the difference between 
longitudinal and transversal structure functions has been 
offered by Siefert & Peinke JgP], who find different cas- 
cade times for longitudinal and transversal spectra. 

As stated above, our results for the longitudinal struc- 
ture functions of the Elsasser variables follow the general- 
ized She-Leveque formula with codimension C — 2 quite 
well. It is difficult to make precise comparisons with ear- 
lier work where often somewhat different cases are con- 
sidered. The perhaps closest comparison is possible with 
the transonic hydromagnetic turbulence simulations of 
Padoan et al. |63 ■ They consider only velocity structure 
functions, but since they have a weak magnetic field, the 
velocity is similar to the Elsasser variables. They find 
that for small Mach numbers the structure function scal- 
ing has codimension C ~ 2, while for supersonic turbu- 
lence they find C « 1. For increasing Mach numbers 
they see a continuous decrease from C = 2 to C = 1. 
This is consistent with our result of C = 1.85 since we 
have a small but finite Mach number, and would there- 
fore expect a somewhat smaller value than C = 2, but 
considerably larger than C — 1. 

In another set of forced turbulence simulations, Cho et 
al 33] find that the velocity structure function has a gain 
codimension C = 2. We emphasize that in Ref. |33 
the structure functions are calculated perpendicular to 
the local magnetic field, and not along the global coordi- 
nate axis. In addition they look at velocity scaling while 
we concentrate on Elsasser variables. Their results can 
therefore not straightforwardly be compared with ours, 
and the similarity is probably accidental. Indeed, it is 
shown in Ref. 66] that structure functions calculated 
perpendicular to the local magnetic field and those calcu- 
lated along the global coordinate axis are not comparable. 

On the other hand, in simulations of decaying magne- 
tohydrodynamic (MHD) turbulence at a resolution up to 
512'^ collocation points, Biskamp & Miiller 17] found the 
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FIG. 19: Longitudinal structure function exponents for the 
Elsasser variables, compared with those for velocity and mag- 
netic field separately for run D. 

codimension of the Elsasser variables to be C = 1. At 
first glance this seems to contradict our results, but the 
fact that they have considered decaying turbulence, and 
that the magnetic field is decaying more slowly than the 
kinetic energy, implies that they have a much stronger 
magnetic field relative to kinetic field than we do. Fur- 
thermore we find (see Fig. I19|) that the magnetic field is 
much more intermittent than the velocity field . It is 
therefore not surprising that they find different intermit- 
tency in their decay simulation than we do in our forced 
simulation. 

Next we look at the longitudinal structure function 
scaling exponents for magnetic field and velocity, as well 
as for the Elsasser variables; see Fig. ^| The velocity 
is known to be generally less intermittent than the mag- 
netic field [g^- While the Elsasser variables follow the 
standard She-Leveque scaling (C = 2) rather well, the 
velocity is less intermittent (C = 3.5) and the magnetic 
field is more intermittent (C ~ 1.7). Obviously, a codi- 
mension larger than the embedding dimension does not 
have a direct geometrical meaning, so the value of C for 
the velocity should not be interpreted as a codimension, 
but just as an indicator for a low degree of intermittency. 

I. Visualizations 

In view of the discussion on the dimensionality of the 
dissipative structures in hydromagnetic turbulence, it is 
desirable to obtain an estimate simply by visual inspec- 
tion. In Fig. 1201 we show, for a small subvolume of 
the entire simulation domain, surfaces of constant 
(Joule dissipation) and (viscous dissipation). Gen- 
erally, the approximate dimensionality of both dissipa- 
tive structures is somewhere between sheets and tubes, 
although appears to be perhaps slightly more sheet- 
like. We conclude that the dimensionality of is con- 
sistent with what one would have expected from the esti- 
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FIG. 20: Contours of (blue or dark-gray) and (green 
or light-gray) showing the dissipative structures of magnetic 
and kinetic energy, respectively. Only a very small subvolume 
(1/16^) of the entire simulation domain of run E is shown. 




FIG. 21: Contours of (blue or dark-gray) and uP' (green or 
light-gray) for the same small subvolume as shown in Fig. 1201 
Note that the areas with strong are more extended than 
those with strong (cf. Fig. 1201 and that the locations of 
high up' and high are near to each other, but otherwise dif- 
ferent in their detailed appearance (see Fig. I2fl^ , even though 



mate C = 1.7 obtained in the preceding subsection. For 
the dissipative structures of kinetic energy, on the other 
hand, the estimate C — 3.5 exceeded the embedding di- 
mension and did therefore not make geometrical sense 
anyway. The qualitative inspection of would have 
suggested a codimension between 1 and 2. We therefore 
conclude that C can only be regarded as a fit parameters 
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FIG. 22: Magnetic field vectors of run E shown at those locations where \B\ > 3-Brms- Note the long but thin arcade-like 
structures extending over almost the full domain. The structures are sheet-like with a thickness comparable to the resistive 
scale. 



and that there is not always an obvious connection with 
the actual appearance of the dissipative structures. 

It should be noted that when one normally talks about 
structures in turbulence one often talks about vortex 
tubes and, in the magnetic case, magnetic flux tubes. 
These are quite distinct from the dissipative structures. 
Vorticity u; = V x u and rate of strain tensor S character- 
ize respectively the antisymmetric and symmetric parts 
of the velocity gradient matrix, and are therefore not ex- 
pected to look similar. On the other hand, (u;^) = (2S^), 
where and angular brackets denote volume averages, and 
Fig.l21lshows that both and exhibit similar length 
scales. The difference is more pronounced in the mag- 
netic case, where J is clearly dominated by smaller scale 



structures while B can exhibit structures of much larger 
scale. This is shown in Fig. 1221 where we visualize mag- 
netic field vectors in the full box at those locations where 
the field exceeds three times the rms value. The strong 
field turns out to be of surprisingly large scale, even 
though this dynamo has no helicity and no large scale 
field in the usual sense; cf. Ref. 5]. The large scale struc- 
tures we find are reminiscent of the ropes discussed in 
Refs [63, 01 , although they seem to be more sheet-like 
|67| . and their thickness is of the order of the resistive 
scale [6^ . 
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nonhelical forcing using the modified forcing function 







FIG. 23: Three-dimensional kinetic and magnetic energy 
spectra, normalized by |«?ma/^i: for runs with and without 
helicity (left and right hand columns, respectively) and for the 
different forcing wavenumbers: kf — 1.5 (top) 2.3 (middle), 
5.1 (bottom). 128^ meshpoints. 



V. CONNECTION WITH HELICAL 
TURBULENCE 



Finally we comment on some differences between heli- 
cal and nonhelical turbulence. In the case of helical forc- 
ing one expects an inverse cascade to smaller wavenum- 
bers, rather than a direct cascade to larger wavenumbers. 
We can now identify two reasons why this has not really 
been seen in early turbulence simulations with helical 
forcing . On the one hand the inverse cascade takes a 
resistive time to develop ^ , and this time tends to be too 
long if magnetic hyperdiffusivity is used . In Ref. 
magnetic hyperdiffusivity was indeed used and the resis- 
tive time was at least a hundred times longer than the 
duration of the runs, so no inverse cascade should be 
expected. But there is another more important reason. 
In order for the inverse cascade to develop, one has to 
have some scale separation, i.e. the magnetic field must 
be allowed to grow on scales larger than the forcing scale 
(which corresponds to the energy carrying scale) of the 
turbulence. This was not the case in the early simulations 
and may explain why the inverse cascade has not been 
seen in Ref. and that those results should therefore 
be closer to the case without helicity. To substantiate 
this, we have carried out simulations with helical and 



fk 



/(nohcl) 
k 



with R," 



(19) 



where jg |;j^g non-helical forcing function of 

Eq. {Tj). In the helical case (a = ±1) we recover the 
forcing function used in Ref. j^, and in the non- helical 
case (cr = 0) this forcing function becomes equivalent to 
that of Eq. 10. 

We show in Fig. the energy spectra of helical and 
nonhelical simulations with forcing at wavenumber kf = 
1.5 (no scale separation), fcf = 2.3 (weak scale separa- 
tion), and fcf = 5 (considerable scale separation). For 
fcf = 1.5 the spectra of the helical and nonhelical sim- 
ulations are indeed quite similar to each other (e.g., no 
inverse cascade and slight super-equipartition at fc > 5). 
For fcf — 5.1, on the other hand, the spectra are quite 
different and there is no inverse cascade in the nonhelical 



VI. CONCLUSIONS 

In the present paper we have studied non-helical 
MHD turbulence without imposed large scale fields. We 
find that, in order to get dynamo action, the mag- 
netic Reynolds number has to exceed the critical value 



Re' 



(crit) 
M 



30. When the dynamo is in the kinematic 
regime (the magnetic field is still too weak to affect the 
velocity) we see a kinetic energy spectrum with the Kol- 
mogorov k~^/^ inertial range, while the magnetic energy 
spectrum shows the expected fc'^/^ Kazantsev slope (in- 
creasing with fc). 

As the dynamo gets saturated we find that there is a 
short inertial range where the magnetic and kinetic en- 
ergy spectra are parallel. This is only seen in the largest 
of our simulations with 1024'^ meshpoints. The magnetic 
energy spectrum exceeds the kinetic energy spectrum by 
a factor of w 2.5, which seems to be more or less the 
asymptotic value as Re grows larger. At first glance one 
is led to believe that the saturated energy spectra exhibit 
a fc~^/^ inertial range, but we argue that this is due to 
a strong bottleneck effect. For one-dimensional spectra 
the bottleneck effect becomes much weaker and they have 
the expected k~^/^ slope. We have demonstrated that 
simulations without magnetic fields show the same bot- 
tleneck effect. Also the second order structure functions 
are consistent with a fc~^/'^ scaling. We therefore con- 
jecture that for larger values of Re one will see a k~^^^ 
subrange also for the three-dimensional MHD spectra, 
although the bottleneck effect will continue to affect at 
least one decade or more in wavenumbers just before the 
dissipative subrange 

At large scales (1 < fc < 5) we see a k^^^ behavior of the 
magnetic field. However, if we force the flow at fcf ^ fci 
(i.e. the energy is being injected at scales much smaller 
than the box size), we find that for fc <C fcf the kinetic and 
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magnetic energy spectra scale almost like fc^, although 
the kinetic energy spectrum is somewhat shallower which 
may be an artifact of the finite size of the computational 
domain. 

Concerning the structure functions we find that our 
simulations are in good agreement with those found by 
Padoan et al. !65j. The Elsasser variables follow the She- 
Leveque scaling with a codimension somewhat less than 
2, which is what one would expect since our simulation 
is weakly compressible. We also find that the magnetic 
field is more intermittent than the velocitVjWhich is qual- 
itatively consistent with earlier findings [lOj. Quantita- 
tively, in terms of structure function exponents, this has 
recently also been found by Cho et al. but it is still 
not known what is the cause of this difference between 
the intermittency of magnetic and kinetic fields. 

In the case of magnetic Prandtl numbers larger than 



unity, there are indications of a range for the mag- 
netic energy spectrum below the viscous cutoff wavenum- 
ber. In order to ensure that this is really the asymptotic 
slope, yet larger simulations are required. 
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